Complementarity and diversity in a soluble model ecosystem 
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Complementarity among species with different traits is one of the basic processes affecting bio- 
diversity, defined as the number of species in the ecosystem. We present here a soluble model 
ecosystem in which the species are characterized by binary traits and their pairwise interactions 
follow a complementarity principle. Manipulation of the species composition, and so the study of 
its effects on the species diversity, is achieved through the introduction of a bias parameter favoring 
one of the traits. Using statistical mechanics tools we find explict expressions for the allowed values 
of the equilibrium species concentrations in terms of the control parameters of the model. 



One of the main interests in the study of model ecosys- 
tems is to elucidate the rules governing the assemblage of 
ecological communities [jjj. In particular, understanding 
the critical role that species composition plays in ecosys- 
tem processes can provide useful guidelines to manage- 
ment efforts on behalf of endangered species (|]. The- 
ory has been prominent in the study of long-term effects 
of species composition on ecosystem functioning, owing 
mainly to the difficulty to keep controlled experimental 
conditions in place for a long period (for a recent remark- 
able exception, see ||). In these studies, N species are 
considered to interact in a given community so that their 
population numbers or biomasses, denoted by x%, . . . , a; at, 
are determined by the dynamical equations 



dxi 



= XiG 
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where G( ) is a species-independent nonlinear function, 
such as Ricker dynamics Q, giving the rate of increase 
of each population and is the interaction coefficient 
measuring the effect of species j on species i (|, ■ As a 
rule, the number of species (diversity) N is taken as the 
control parameter or independent variable and the total 
biomass (productivity) Xi, viewed as a measure of the 
ecosystem stability, as the dynamical variable (|, |§[ 0- 
However, especially in long-term studies of large ecosys- 
tems this approach is faulty since it does not take into 
account the dynamics of diversity and the species poten- 
tial for adaptation to environmental fluctuations 0. 

In this contribution we consider an alternative ap- 
proach in which the diversity is given by the number of 
surviving species in the ecosystem at equilibrium, being 
thus a dynamical variable that, ultimately, depends on 
the nature of the species composing the ecosystem. In 
particular, we assume that the fraction or concentration 
of individuals of species i in the ecosystem, Xi € [0, oo), 
is determined by the nonlinear system of equations, so- 
called replicator equations 0, 



dt 



(2) 



where T% = Ylj Jij x j can be identified with the fitness of 
species i and the term <f> = ^\ Xi Ti ensures that 



= N 



(3) 



for all times. This constraint enforces an effective com- 
petition among the species in the ecosystem. We note 
that by dropping the <f> term in Eq. (^|) we recover the 
classical equation of Lotka- Volterra Q . 

In the case of symmetric interactions = Jji the 
asymptotic regime of Eq. (^|) is simply characterized: 
the dynamics maximizes the Lyapunov function or fitness 
functional 



JijXiXj 



(4) 



and so it can be shown that the only stationary states are 
fixed points || . In this case Jy < corresponds to pairs 
of cooperating species whereas Jy > to pairs of com- 
peting species. The equilibrium regime as well as some 
aspects of the dynamics can become nontrivial, however, 
in the case that the functional T has many local maxima. 
This occurs when the coupling strengths are quenched 
random variables taking on positive and negative values, 
as in the model of random replicators put forward by 
Diederich and Opper Jl(J (for further development see 
p| |l^|). The assumption of random interactions is a 
form of taking into account our lack of knowledge of how 
the species actually interact. Moreover, it allows the use 
of tools of the statistical mechanics of disordered systems 
to fully characterize the equilibrium states, thus making 
feasible the study of large ecosystems. In contrast, to 
keep numerical accuracy under control, the traditional 
approach based on the numerical solution of Eq. ([!]) is 
restricted to ecosystems composed of typically N = 10 
species. We note that in the random replicator frame- 
work the productivity J2i x i i s constant, while the diver- 
sity varies since a fraction of the N species may go extinct 
due to outcompetition. This phenomenon becomes ap- 
preciable, and hence passive of quantitative analysis, for 
large N only. 
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Up to now studies of the random replicator model 
have considered the strengths of the interactions be- 
tween species as independent, Gaussian distributed ran- 
dom variables |L0[ O, [l2|. In this contribution we go 
beyond this initial stage by introducing some underlying, 
non-random structure in the interspecies interactions. In 
the spirit of models for molecular recognition Jl3|, |l4j we 
assume that each species is characterized by a set of p 
traits, £f,/i = l,...,p, and that the resulting interac- 
tions between pairs of species depend on these traits ac- 
cording to a complementarity principle. Specifically, we 
assume that the traits £f are quenched, independent ran- 
dom variables that can take on the values +1 and —1 with 
probabilities (1 + a)/2 and (1 — a)/2, respectively. Here 
a e [0, 1] is a bias towards the trait +1, corresponding, 
e.g., to an economically favored feature of the species. In 
this sense, the parameter a may be thought of as a mea- 
sure of the human impact on the species composition of 
the ecosystem. In addition, we assume that the coupling 
strength between species i and j is given by the Hebb 
rule 



p 



i + 3 



(5) 



that was extensively studied in the eighties within the 
neural networks context P, |J|. Clearly the lar ger the 
number of complementary traits (i.e., £,f£^ = — 1), the 
more cooperative the pair of species. Complementarity 
among species with different traits has been suggested 
as one of the major mechanisms involved in biodiversity. 
In that context, the term complementarity subsumes all 
local deterministic processes which increase the perfor- 
mance of communities above that expected from the per- 
formance of individual species grown alone ||. If com- 
munity performance is measured in terms of the fitness 
functional (||) then the Hebb coupling between species 
(|J) is clearly well suited to model such processes. 

As the final ingredient to define the model, we need 
to specify the self-interactions J a. Though for finite N 
the constraint (|3|) prevents the unbounded growth of any 
single species, it becomes inefficacious in the thermody- 
namic limit N — > oo and so an additional mechanism to 
limit growth becomes necessary. The usual procedure is 
to introduce a competition term between individuals of 
a same species by setting Ju = u > for all i so that 
u can be viewed as a global cooperation pressure. In 
practice, a positive self-interaction is essential to ensure 
convergence of the numerical methods used to solve the 
replicator equations @). 

The maxima of the functional given in Eq. (^) can 
easily be obtained within the statistical mechanics frame- 
work in the limit of infinite N but finite p (see, e.g., [[l6f ). 
We begin by defining the free-energy density / as 



lim — In Z 
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where 

Z = 



POO 

Jo , i 

fj, \ i / i 



x exp 
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is the partition function and j3 = 1/T is the inverse tem- 
perature. Taking the limit T — > in Eq. (6) ensures 
that only the states that maximize T will contribute to 
Z. After some standard algebric manipulations |l6[ , we 
find 



Pf = 



_ m- , 1 [7 



(Q + m-0' 



In crfc 



Q + m-£ 

2y^ 



(8) 



with the notation z ■ y = z^y^ and where we 
have invoked the self-averaging property of the sums 
X)jS(£f)> where g is any continuous function, 
to replace them by the species- independent averages 
(#(£ M ))- Here (...) stands for the average taken with 
the probability distribution 



p(o=n 



1 + a 



s(e-i) + 



1 



(9) 



The saddle-point parameters Q and are given by the 
solutions of the equations df/dQ — and df/dm^ = 
0. Only has a relevant physical meaning, namely, it 
is the average overlap between the equilibrium solutions 
and the trait £f 



1, 



,N, 



m' 



(10) 



where (. . .)t stands for a thermal average taken with the 
probability distribution 

W({z 4 }) = U{N - Xi) exp \pF({xi})] . (11) 



The next step is to take the zero-temperature limit /3 — > 
oo. Nonzero solutions for Q and m M are found in the 
regime Q + m • £ < only, where the p + 1 saddle-point 
equations take on a particularly simple form, whose only 
solution is the symmetric one m M = m V/x. Explicitly, we 
find 



2ua 



2u+l- a 2 



and 



= -2u 1 



pa 



2u + 1 - a? 



(12) 



(13) 
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with the notation Q = Q/(3 and rh = m/(3. 

The characterization of the ecosystem through the 
global parameter m is not very illuminating, and a better 
understanding can be achieved by looking directly at the 
values the species concentrations Xi can take on. This 
can be done, for instance, by calculating explicitly the 
cumulative distribution that the concentration of a given 
species, say Xk, assumes a value smaller than x, defined 
by 

C k (x)=\im / TJdxj e(x-x k )W({xi}) (14) 

where &(x) = 1 if x > and otherwise, and W ({xt}) 
is given by Eq. ([ll]). Since all species concentrations are 
equivalent we can write Ck (x) — C (x) Vfc and evaluate 
Eq. ( |l4| ) by adding the field term h Q (x — Xi) to Eq. 
(^). Taking the derivatives of the resulting free-energy 
with respect to h and then the limit h — > yield 



C(X) = J2 Wn& (X ~ Cn) , 



n=0 



where 



and 
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Hence the species concentrations can take on the val- 
ues x = Cn (provided that Cn > 0) only, and the frac- 
tion of species with concentration Cn is given by W n . 
Here the integer variable n yields the total number of 
traits +1 assigned to a given species. The fact that 
ii is the only species feature that determines the equi- 
librium concentration is consequence of the symmetric 



solution 



= m V/i of the saddle-point equations. 



Since Co > Ci > • • • > Cp the economically more rele- 
vant species, i.e., those characterized by the set of traits 
C^ = +1 V/i, are the first ones to go extinct when, say, the 
cooperation pressure u decreases. The condition for the 
coexistence of all species in the ecosystem is then Cp > 
that reduces to 



p < 1 



1 



2u 



a (1 — a) 



(18) 



We note that the species characterized by n < p/2 never 
die out. A comparison of these analytical results with 
the numerical solution of Eq. (^) for N = 2000 and a 
single instance of the interaction matrix is presented in 
Fig. ^] where we show the cumulative distribution C{x) 
for p = 5. For the sake of illustration, the corresponding 
probability distribution P(x) = dC/dx is also shown in 
the inset. 




FIG. 1: Cumulative distribution of species for p = 5, u = 0.1 
and a — 0.5. The dashed curve is the prediction of Eq. ( |li| ) 
and the solid curve is the result of the numerical solution 
of the replicator equations. The inset shows the histogram 
associated to the probability distribution P(x). 



The situation of maximum cooperation occurs for a — 
since Cn = 1 Vn and any particular set of p traits con- 
tributes with the same fraction, 2~ p , to the final ecosys- 
tem composition. However, this symmetry is broken for 
a > 0. In fact, as the result of the excess of traits 
C 1 " 1 = +1, many different species carry the same or nearly 
the same set of traits and so their interactions have a 
strongly competitive character, i.e. ~ p/N > 0, ex- 
plaining thus the depletion of the concentration of species 
with n > p/2. On the other hand, the advantage of 
species characterized by an excess of traits £ M = — 1 is 
twofold. First, their relative rarenesses imply that there 
are few competing species and, second, they have abun- 
dance of cooperative partners among the species with 
n > p/2. These two factors explain the dominance, in 
the sense of occuring at larger concentrations and hence 
of possessing greater stability, of the species character- 
ized by n < p/2. A particularly extreme situation is 
observed in the limits M-t0 and a — > 1, when all species 
with n > p/2 die out. These findings can be given an ob- 
vious ecological interpretation, namely, species that are 
too similar and, consequently, compete for essentially the 
same resources are proner to extinction. However, the 
situation is not so simple in the case of a nonzero u. 
For instance, in the limit a — > 1 the explosive growth of 
species with n < p/2 can be controlled even by a van- 
ishingly small self- interaction parameter, resulting in the 
coexistence of all species [see Eq. (|l8|)]. Fig. |2| illus- 
trates the complex dependence of the species diversity 
D = 1 — C(0) on the control parameters a and p for 
u = 0.1. The results for u = are very similar, except in 
the close neighborhood of a = 1 where D — > 0, implying 
that there are a finite number of surviving species only. 
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correlation between diversity and productivity in exper- 
iments on randomly assembled communities, whereas in 
nature the most productive ecosystems are those charac- 
terized by low species diversity B. These antagonistic 
conclusions may well be due to the attempt to describe 
productivity solely in terms of diversity, while a complete 
description would require the knowledge of other, prob- 
ably uncontrolled, quantities. 
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FIG. 2: Fraction of surviving species or diversity D as function 
of the bias parameter a for u — 0.1 and p = 5 (O)i 10 (v); 
and 15 (A). The solid curves are the theoretical predictions 
and the symbols are the results of the numerical solution of 
the replicator equations. 

These results indicate that, when complementarity is 
the sole mechanism determining the species interactions, 
the assemblage of synthetic ecosystems aiming at the ex- 
ploitation of some particular traits (as in monocultures, 
for instance) may, in the long term, be disastrous to the 
economically relevant species. However, if there is in 
addition some external global pressure for cooperation, 
then the best strategy for the long-term survival of those 
species is to guarantee that they are massively present in 
the initial assemblage. Even so, the more stable species 
are always those of less economic value. 

Although diversity has been manipulated as an in- 
dependent variable in most of the numerical studies of 
model ecosystems as well as in small scale experiments, 
it is becoming evident that, in large scales, species di- 
versity itself is a dynamical variable that adjusts freely 
to changes in environmental conditions ||. The under- 
standing of more complex ecological systems calls then 
for a more holistic approach. Though relying on the un- 
realistic, technical assumption of symmetry of the inter- 
actions between species, equilibrium statistical mechan- 
ics provides an useful (if not the only) analytic frame- 
work to tackle this difficult issue. The elusive character 
of diversity is apparent from our results (see, e.g., Fig. 
^) since, despite the simplicity of our model ecosystem, 
its complex dependence on the control parameters of the 
model precludes a simple description in terms of a sin- 
gle parameter. This situation is reminiscent of a recent 
impass in ecology, triggered by the finding of a positive 
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